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We revisit the problem of a triad of resonantly interacting nonlinear waves driven by an external 
force applied to the unstable mode of the triad. The equations are Hamiltonian, and can be reduced 
to a dynamical system for 5 real variables with 2 conservation laws. If the Hamiltonian, H, is zero 
we reduce this dynamical system to the motion of a particle in a one-dimensional time-independent 
potential and prove that the system is integrable. Explicit solutions are obtained for some particular 
initial conditions. When explicit solution is not possible we present a novel numerical/analytical 
method for approximating the dynamics. Furthermore we show analytically that when H = the 
motion is generically bounded. That is to say the waves in the forced triad are bounded in amplitude 
for all times for any initial condition with the single exception of one special choice of initial condition 
for which the forcing is in phase with the nonlinear oscillation of the triad. This means that the 
energy in the forced triad generically remains finite for all time despite the fact that there is no 
dissipation in the system. We provide a detailed characterisation of the dependence of the period 
and maximum energy of the system on the conserved quantities and forcing intensity. When H ^ 
we reduce the problem to the motion of a particle in a one-dimensional time-periodic potential. 
Poincare sections of this system provide strong evidence that the motion remains bounded when 
H 7^ and is typically quasi-periodic although periodic orbits can certainly be found. Throughout 
our analyses, the phases of the modes in the triad play a crucial role in understanding the dynamics. 

PACS numbers: 05.45.-a,02.30.Ik,92.10.hf 



I. INTRODUCTION 

Resonant triads are the basic building blocks for under- 
standing mode coupling in systems of weakly interacting 
dispersive waves in which the leading order nonlinearity 
of the underlying wave equation is quadratic in the wave 
amplitude, "0. Physical examples of such systems include 
capillary waves on fluid interfaces [TJ, Rossby waves in 
geophysical fluid dynamics [3] , mode coupling in non- 
linear optics [1], drift waves in magnetized plasmas [3] 
and internal waves in stratified fluids [7] . For simplic- 
ity, let us suppose we can represent such wave fields in 
terms of Fourier harmonics. Each linear mode can be la- 
belled by its wave- vector, k <E M. d (where d is the physical 
dimension of the system), and frequency, Wk- A resonant 
triad is a triplet of modes, ki, k2 and k3 which satisfy 
the resonance conditions: 



ki = k 2 + k 3 



(LI) 



In the limit of weakly nonlinear waves, energy is only 
exchanged efficiently between modes which are members 
of resonant triads. For an idealised wave system in an 
unbounded domain, the k's are continuous variables. A 
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consequence of this fact is that a given k is typically a 
member of many resonant triads. The resonant triads 
therefore form a network in Fourier space which provides 
a conduit for energy transfer throughout the full set of 
Fourier modes. This idea of a network of resonant modes 
carrying energy throughout the system is the principal 
basis for the theory of wave turbulence [5] . Wave turbu- 
lence describes the statistical dynamics of an ensemble of 
weakly interacting dispersive waves in which many modes 
are excited by external forcing. In the usual wave turbu- 
lence scenario, which is relevant in many applications, the 
forcing directly supplies energy only to a subset of these 
modes, typically those corresponding to large scales. En- 
ergy is then transferred to other modes, in particular to 
small scale modes, via resonant interactions in a process 
referred to as an energy cascade. The theoretical de- 
scription of such weakly nonlinear cascades is very well 
developed (see [5] for a recent review). 

For wave systems in bounded domains several inter- 
esting complications arise. The k's are now discrete vari- 
ables. In the simplest case, of a multiply periodic domain 
of size L, the k's are of the form k^ Ak where Afc = 2 ir/L 
and krf £ Z d is a vector of integers. Finding resonant 
triads in discrete systems is therefore equivalent to find- 
ing integer-valued solutions of the resonance conditions, 
Eq. (1.1). In many applications, the dispersion relation, 



Wkj is a power law, Wk = k a , often with a fractional ex- 
ponent, a. For such dispersion relations, integer valued 
solutions of the resonance conditions, if they exist at all, 
are quite rare. Capillary waves (a = 3/2, d — 2) are 
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a well known example for which no discrete solutions of 
the resonance conditions are possible [10] • The sparsity of 
discrete solutions of the resonance conditions means that 
any particular k in a bounded wave system is typically 
only a member of relatively few resonant triads. Triads 
which share a common mode can be grouped together 
to form resonant clusters. While many such clusters of 
varying size coexist in the Fourier space, no one cluster 
spans the entire space and the resonant energy trans- 
fer network fragments into disconnected parts. A full 
enumeration of the resonant clusters for Rossby waves is 
studied in [TT] . While very large clusters are typically few 
in number, there is usually a large number of small clus- 
ters. In particular, isolated triads are quite common. An 
isolated triad is a triad which shares no common mode 
with any other triads. If the energy of the system is ini- 
tially concentrated in the modes of a small cluster, the 
simplest case being an isolated triad, then no resonant 
energy cascade can take place. This gives rise to the 
phenomenon of "frozen turbulence" [12j in which energy 
remains localised in relatively few modes rather than cas- 
cading throughout the system. In practice, the level of 
nonlinearity is always finite leading to nonlinear broad- 
ening or detuning of resonances. This finite amplitude 
effect then allows energy to be exchanged between triads 
which are not exactly resonant so that a cascade may 
proceed even in the absence of a large resonant cluster 
spanning the system [13] . Such quasi- resonant interac- 
tions are less efficient at transferring energy in the sense 
that the timescale for energy transfer is longer compared 
with resonant interactions. For this reason, the exactly 
resonant clusters, in particular the isolated triads, remain 
the building blocks for understanding energy transfer in 
the system when the level of nonlinearity is finite but 
small. 

After suitable rescalings of the wave amplitudes, the 
equations for an isolated triad can always be brought to 
the canonical form: 





— ZB^Ba, 




B 2 


= ZB*B 3 , 


(1.2) 


B, 


= —ZBiB 2 , 





where B\, B 2 and B% are the rescaled complex ampli- 
tudes of the modes constituting the triad. This is a 
Hamiltonian system with Hamiltonian 



H = 2Zlm{B 1 B 2 B^}. 

With this Hamiltonian, Eqs. ( |I.2[ ) are equivalent 
Hamilton's equations: 

8H 



to 



iB, 



dB* 



J = 1,2,3 



along with their complex conjugates. H is obviously con- 
served by these dynamics. 

Eqs. (1.2 1 have been extensively studied for many years 



since early work in nonlinear optics [4], plasmas, atmo- 
spheric dynamics |14j and electronics led to an appreci- 
ation of the importance of mode coupling in nonlinear 



systems. They describe the nonlinear saturation of the 
so-called decay instability which manifests itself in a vari- 
ety of wave systems with quadratic nonlinearity including 
capillary waves [Tj , Rossby waves [2J [3] and drift waves 
[5]. The decay instability, in its linear stage, describes a 
process whereby modes B\ and B 2l which initially con- 
tain little energy, grow in amplitude at the expense of 
mode -B3 which initially contains most of the energy. For 
this reason, B3 is referred to as the unstable mode of the 
triad. More recently, Eqs. (1.2) have been studied in a 



mechanical context to describe the dynamics of a spring 
pendulum (T51I16) . Moreover there has been considerable 
recent interest in using these equations as models of reso- 
nant Rossby wave interactions to attempt to explain the 
observed periods of various intra-seasonal oscillations in 
the Earth's atmosphere [I?] . 



It has been known for a long time that Eqs. (1.2) form 
an integrable system (for a discussion see [18J ) . Analytic 
formulae for the amplitudes of the modes in terms of 
elliptic functions have been known at least since the early 
1960's |4]. Recently, corresponding explicit formulae for 
the phases, as encoded by the dynamics of the so-called 
dynamical phase, were obtained in [Tj5]. This work means 
that everything is effectively known analytically about 
the dynamics of an isolated triad. 

Given the role played by resonant triads in the theory 
of wave turbulence discussed above, it is somewhat sur- 
prising that almost nothing is known about the dynamics 
of an isolated triad in the presence of forcing. We aim to 
address this in the present paper. We consider the sim- 
plest possible forcing: an additive forcing of the unstable 
mode. We force the unstable mode because this will al- 
low efficient transfer of energy into the triad since the 
unstable modes tends to transfer energy to the other two 
members of the triad whereas energy in the stable modes 
tends to remain there. Besides its ubiquity in the study 
of wave turbulence at both the numerical and analytic 
level, the addition of forcing to the equations of motion 
has direct physical meaning in certain applications. In 
the case of Rossby waves, for example, an additive forcing 
term in Eqs. (1.2 1 describes orographic forcing of Rossby 
waves by an idealised periodic topography [5D]. For the 
swinging spring it would correspond to a periodic oscil- 
lation of the point of support. A numerical study |21j of 
the forced triad equations suggested that the amplitudes 
of the modes remain bounded in the presence of forcing, 
something which is not obvious a priori. 

In this paper we provide a detailed study of the prop- 
erties of an isolated triad in the presence of forcing. The 
paper is laid out as follows. We begin in Sec|TT]by writing 
the equations of the forced triad in the general case and 
studying their conservation laws. We then specialise in 
Sec. [HI] to studying the case H = 0. We prove that the 
system is generically bounded when H — meaning that 
it is bounded for any choice of initial condition with the 
exception of a single special initial condition for which 
the system is unbounded (see Sec. HID 2). We obtain 



an explicit solution in Sec III C of the forced triad evolu 
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tion equations for a particular choice of initial conditions. 
An explicit solution of the evolution equations in terms of 
elliptic functions does not seem possible in general, even 
in the integrable case H = 0. Instead, in Sec. |IV| we 
provide formulae for the period and maximum amplitude 
of the solution which are valid in general. We perform 
a parametric study of the period and maximum ampli- 
tude of the triad as key parameters are varied. Then in 
Sec.|VJ we introduce a novel and accurate scheme for ap- 
proximating the H = dynamics and compare it against 
numerical solutions of the original system for a variety 
of initial conditions. Finally, in Sec. |VI[ we turn our 
attention to the case H ^ 0. We find that the dynamics 
can be reduced to the one-dimensional motion of a par- 
ticle in a time-periodic potential. Numerical Poincare 
sections of the resulting dynamical system suggest that 
the dynamics remain bounded when H ^ and the mo- 
tion is generally quasi-periodic although periodic orbits 
can also be identified. The article then ends with a brief 
summary and conclusions. 



Under this new representation, the equations of motion 



(II.l I become 



f Ci 

c 2 
c 3 

01 

02 
03 



ZC 2 C 3 cos ip, 
ZC\C 3 COS if, 

—ZC\C 2 cos ip + F sin ip 3 , 

ryC 2 C 3 . 

-Z——s\mp, 
-Z——smtp, 

L-2 

7 CiC 2 . F 
03 03 



(II.2) 



where we have introduced the dynamical phase ip = ip.\ ■ 

<P2 - ^3- 



B. Conservation Laws of General Forced Triad and 
reductions 



II. GENERAL FORCED TRIAD 

Let Bi, B 2 and B 3 be complex functions of time. We 
assume they satisfy a system of evolution equations de- 
rived from the real Hamiltonian 

H = 2ZIm{B 1 B 2 B 3 } - 2Flm{iB 3 }, 

where F and Z are constant real parameters. The canon- 
ical equations of motion for this system are given by 



r)ff 



along with their complex conjugates. The equations of 
motion are, therefore, 



The unforced triad (i.e., the case F — 0) is known 
to possess two independent conservation laws that are 
quadratic in the amplitudes Cj , and one conservation law 
which is cubic in the amplitudes. In contrast, the forced 
triad possesses one quadratic conservation law and one 
cubic conservation law. These are 

J = C x — C 2 , 

H = 2C 3 [ZdC 2 sin (p - F cos tp 3 } . (II.3) 

Here, H is nothing but the Hamiltonian of the original 
system (II.l), written in terms of the new variables. 



Another important property of the unforced triad is 
that all the phases tp±, ip 2 and ^3 were 'slave' variables, 
i.e., they were obtainable by quadratures once the solu- 
tions for Ci, C 2l C3 and tp were known. In contrast, in 
the forced triad only ifi and (p 2 are 'slave' variables, and 
the five real variables C\, C 2l C3, ip and ^3 form the 
reduced system 



Bi — ZB^B 3 , 
B 2 = ZB{B 3 , 
B 3 = -ZBxB 2 +iF, 



(HI) 



together with their complex conjugates. In this setting, 
the term 'iF ' is originated from a forcing of the so-called 
'unstable' mode represented by B 3 . 



C\ = ZC 2 C 3 cos tp, 
C 2 = ZC1C3 cos<^, 
C3 = — ZC±C 2 cos ip + F sin (£3, 
H 



<P3 = 



2Cf 



-ZC\C 2 C 3 sinyj 



1 

ci 



1 



H 

2Cl> 



(II.4) 



A. Amplitude-Phase Representation 

We write Bj = Cje^ , j = 1,2,3, with Cj € K and 
ipj € K to obtain 



Bj = (Cj+iCjV 3 )^. 



where H is defined as in equation (II. 3 ). Numerical simu 



lations of system (II.4) can be interpreted as either using 



H as a constant or by explicitly using H as defined in 
equation (II.3). The accuracy might be sensitive to this 
choice. 



System (II. 4 1 along with the two conservation laws, 



J and H , is effectively a three-dimensional system and 
as such it might not be integrable. We notice, however, 
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that the system is volume-preserving, with Jacobi last 
multiplier 



C1C2, 



if H = constant, 



CiC 2 C 3 , if H is defined as in Eq.(II.3) 



Remark: a Jacobi last multiplier, also known as a stan- 
dard Liouville volume density, is a scalar function p(x) of 
the dependent variables x a , a = 1, . . . , N, defined by 
the equation V • (p V(x)) = , equation to be valid for all 
values of x, where V is the right-hand-side of the evolu- 
tion equations x = V(x). See [221 - E4"] . [25] and references 
therein. 



C. Integrability of the case H = 



We consider initial conditions for system (11.4) such 
that H = 0. We obtain 

Ci = ZC2C3COSV2, 

C2 = ZC\C 3 COS If , 

C 3 = —ZC\Ci cos ip + F sin (^3, 
^3 = 0, 

1 1 



(II.5) 



-ZC1C2C3 simp 



cl 



c\ 



We notice that the case C3 = is trivially integrable. As- 
suming that C3 ^ 0, we obtain H = ZC1C2 simp = 



F cos 933. But, according to the dynamics (|II.5|, tp 3 is 
constant, so we derive the new conservation 



C1C2 sin ip 



F cos ip 3 



(II.6) 



This new constant, together with J = Cf — Cf, helps 
us reduce the system to an effectively two-dimensional 
volume-preserving and therefore integrable system: 



' Ci = ZC2C3 cos ip, 
C2 — ZC1C3 cos cp, 
C3 = —ZC1C2 cos ip + F sinipz, 
1 1 



(II.7) 



-FC3COS p> 3 



ci 



with </?3 being constant, together with the conservation 
laws 



C1C2 sin</3 



F cos tp3 



(II.8) 
(II.9) 



and Jacobi last multiplier given by 

Ci C2 , if i 71 cos (p 3 = ZC\ C2 sin ip 

P 



1. 



is used in Eq.(II.7l, 

if F cos tps = constant 



is used in Eq.(II.7l 



III. GENERIC BOUNDEDNESS OF THE 
INTEGRABLE CASE H = 

The goal within this section is to establish generic 
boundedness of the solutions to the integrable system 
(II.7HH.9). By boundedness we mean that the 'energy' 
jUJT) of the system is bounded for all times. The energy 
is a positive-definite quadratic function of the amplitudes 



E(t) 



C 1 (t) 2 + C 2 (tf 



C 3 (t) s 



(III.1) 



This quantity is obtained from the kinetic energy of the 
original wave system and should not be confused with 
the Hamiltonian H, which is a cubic function of the am- 
plitudes and equal to zero in this case. 

The time derivative of the energy is obtained using 
equations ( II.7[ ): 



E(t) = C 1 (t)C 1 (t) + C2(t)C 2 (t) + 2C 3 (t)C 3 (t) 

= 2F sm(p 3 C 3 (t) . (III.2) 

This result suggests that we use a 'local rescaled time' 
variable r(t) satisfying 



dr 
dt 



C 3 (t), r(0)=0. 



(III.3) 



This time transformation is well defined as long as C 3 {t) 
is not zero. Integrating equation (III. 2) we obtain 



E{t) = E(0) + 2tF sin ^ 3 



(III.4) 



where r = r(t). We see that proving boundedness of 
energy is equivalent to proving boundedness of r(t). 



A. Establishing boundedness for J 7^ 

The conservation law J = Cf — Cf naturally suggests 
a parameterisation of the form: 

/ Ci(t) = \/jsgn(Ci(0))coshw(r), 
1 C 2 (t) = \/jsinhw(r), 

in the case J > and 

J Ci(t) = v^J sinh uj{t), 

\ C 2 {t) = \/=Jsgn(C 2 (0))cosh W (r), 

when J < 0. Substituting this parameterised solution 



into equations (II. 7 1 and using equations (II.8I, (II.9) and 



(III.3I, we deduce that the function lo{t) must satisfy, 
of whether 

cosh 2w(t) 



regardless of whether J is positive or not, 

2 



i?cosh (2Zt + A) 
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where R > and A g R are constants, related to the 
initial conditions as follows: 



R = 




F 2 COS 2 (fi3 



R sinhA = Ci(0)C 2 (0) cos^(O). 



(III.5) 
(III.6) 



Substituting the solution for w(r) into the parame- 
terised forms for C\(t) and C2(t) gives the result valid 
for J ^ : 



Ci (t) 2 = R cosh (2Zr + A) + 
C 2 (t) 2 = R cosh (2Zr + A) - 



J 

2" 

J 

2' 



(III.7) 



The solution for the remaining amplitude C 3 (i) is found 
by using equations ( III. 1 ) and (III.4). We get 

C 3 {t) 2 = E(0) + 2Ft sin ^ 3 - i?cosh (2Zt + A0IJI.8) 



and using equation (III. 3 1 we derive the following equa- 
tion for r(t): 



dr 



i?cosh(2Zr + A)-2FTsin(^3 = -E(Q). (III.9) 



This is precisely the equation for a one-dimensional par- 
ticle with position r(t), total energy -E(O), kinetic energy 

dr" " 

d/ 



and potential energy 
V(t) = i? cosh (2Zr + A) - 2Fr sin ip 3 . 



(III. 10) 



From the fact that 72 > Q it follows that this poten- 
tial is in fact attractive: V"(t) > for all r G R, and 
lim T _>±oo ^( r ) = oo. It follows that the particle must 
oscillate about a unique minimum situated at the point 



A 



1 



2Z 2Z 



sinh 



F sin (fi3 
ZR 



with the motion of r(t) being bounded between the two 
turning points T m - m < r max defined by the solutions of 

V(T min ) = ^(w) = E(0), 



Since r(t) is bounded, we conclude from equations (III. 7) 
and ( |IIL8 ) that the amplitudes of each mode must also 
be bounded. 

Evolution of the dynamical phase ip(t). Equation 
(II. 8 1 could at first sight be used to solve for the dynam- 



ical phase ip(i), however this is not practical because the 
dynamical phase generically oscillates about ir/2 (mod- 
ulo mr) and therefore the inverse sine is multivalued. It 
is better to integrate directly equation ( II.7[ ) for (p using 
the auxiliary variable r{t). The result is: if the initial 
condition <f(Q) is in the open interval (rt7r, (n + for 
some n € Z, then ip(t) remains in that interval and 



ip(t) = cot 



_! Zi?sinh(2Zr + A) 

F COS f 3 



(III.ll) 



where the branch is properly selected. In the special case 
when the constant ip 3 is equal to (to + l/2)7r , m G Z, 
then if is also constant and equal to run , m 6 Z . 



B. Formulae for the turning points in terms of new 
special functions 



The turning points r, 



mm / max 



of the evolution equation 



(III. 9 1 are defined by the solutions of the equation V(t) 



R cosh (2Zt + A) - 2Ft sin p 3 = E(0) . 



(111.12) 



For physically admissible values of the parameters Z, F 
and initial conditions R, f 3 ,E(0), this equation has ei- 
ther two real solutions or only one real solution. When 
there is only one real solution, the physical system is 
"dynamically frozen" at r = constant. In order to study 
a less trivial nonlinear dynamics, it is better to restrict 
the parameter space so that there are two real solutions, 
denoted by r min ,T max , with r min < r max . 

Let us study the structure of these solutions as func- 
tions of the parameters. Noti ce that if we replaced the 
hyperbolic cosine in equation ( III.12[ ) by an exponential, 
the solution for r min / max would be given in terms of two 
branches of the product log function Wk(x), a function 
of one variable. In contrast, in our case such reduction 
does not get too far: the solution is given in terms of a 
function of two variables. Defining the function £,(x,y) 
by the implicit equation 



cosh £ = x £ + y , 



(111.13) 



we look for real solutions (eR. Graphically, we want to 
find the intersections between two graphs, cosh£ and 
y, as function of x and y. We see that tangent intersection 
occurs if in addition to the above equation, sinh£ = x 
holds. Therefore, a relation between the arguments is 
found in the tangent case: y = \/l + x 2 — x sinh -1 x . 
Therefore, the condition (of physical origin) to have two 
solutions is a restriction on the parameter space: 



V > 



x sinh x , x G 



(III. 14) 



With this condition, there will be two real solutions of 



equation (III. 13), denoted 

With these functions defined, the turning points 
fmin, T max can be obtained via proper rescalings: 



2Zt„ 



mm / max 



A = £ fc 



Fsin^ 3 E(0) F sin c^ 3 A 



R Z R R Z 

(111.15) 

where k = or — 1. In practice, we resort to a numerical 
algorithm to compute these solutions. 

Condition (III. 14 1 is automatically satisfied for any 



choice of real initial amplitudes and phases. 



() 



C. An explicit example of bounded motion: 

sin (£3 — 0, J ^ 

Continuing with the case H = 0, consider now a choice 
of initial conditions such that ip3 £ {0, ±7r} and J ^ 0. 
The differential equation for r(t), equation (III. 9), now 
reads 



dr 
dt 



= E(0) - R cosh (2Zt + A) 



(111.16) 



The particle interpretation is simple, given the fact that 
E(0) > R by virtue of the condition H = 0. It turns out 
that the potential V(t), apart from a shift, is a symmetric 
well and therefore the motion should be bounded. This 
equation can be explicitly integrated to give an analyti- 
cal solution for r(t) in terms of Jacobi elliptic functions, 
with modulus m = ^|°| + ^ € (0, 1) . The procedure is 
standard but tedious, and is omitted here. The result can 



be checked directly by substitution into equation (III. 16 1. 
We obtain the explicit solution 



2Zr{t)+A = In 



[d„(: 



■ £+u|m)- 



2 JC(m) 



1 





where K(m) is the complete elliptic integral of the first 
kind, the period T of the oscillations is defined by 



T = 



2K(m) 



Zy/E(0)+R' 



and the shift u is defined in terms of the incomplete el- 
liptic integral of the first kind by 



u=f( g'( > 



To see that the above solution for r(t) is indeed oscil- 
latory, we compute its time derivative: 



dr r —, (2 Kim) 

— = ^E(0)-Rsn(^ r ^t + u 



With this result, we readily obtain the squares of the 
amplitudes C 1 (t),C 2 (t),C 3 (t) : 



, 2K(m) 
^ivJ = e(0)+£ - (e(o)-r) sn" 1 

C 2 (t) 2 = E(o)-i - (E(O)-R) sn 2 



T 

2K(m) 



T 



t + u 
t + u 



C 3 {t) 2 = (E(o)-R) sn 2 f A ' t + ■ 



We plot these squares in figure[T] top panel. Finally, from 
equation ( III. 1 1 ) we obtain the dynamical phase: 



¥>(*) 



cot 



Fl f R f /2 cndnf 2 -^t + u 
F (1— m) cos ^3 V T 



mj 
(111.17) 



FIG. 1. Solutions for the square amplitudes (top panel) and 
dynamical phase (bottom panel) in the explicitly integrable 
case H = and tp3 — 0. Choices for parameters and initial 
conditions: F = 1, Z = 1, E(0) = 3, J = 1, and ip(0) is 



obtained from equation (III. 17 1 



where we have defined cn dn(x|m) = cn(a;|m) dn(a:|m). 
Depending on the initial value f(0), the dynamical phase 
remains bounded and oscillates about the value (n + 
l/2)7r, for some n G Z. The corresponding plot of the 
solution for the dynamical phase is shown in figure [l] 
bottom panel. 



D. The limit J — > 0: "almost always" bounded 

Always in the case H — 0, we have established bound- 
edness of the motion for the case J 7^ in Section III A 



Using the equations from that Section, it is possible to 
establish boundedness in the case J = 0, cos ip 3 ^ 0, by 
taking the limit of the J ^ equations as J — > 0. For 
example, the auxiliary function w(r) has a singularity in 
the limit: 



u(t) ~ ^log 



4|Fcos ipz\ 
\JZ\ 



cosh(2Zr + A) 



but the physical quantities such as amplitudes and phases 
behave well: in particular, equations (III.7I and (III. 11) 
hold, with the only changes J 



and R = lFc ° z s ^ . 
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What about the case J = and cos ^3 = 0? To get 
there by a limiting procedure of the J / equations, 
we need to be more careful. The condition J = im- 
plies Ci(i) 2 = C2(t) 2 , Vt > 0. For simplicity, we take 
initial conditions Ci(0) = C 2 (0) ^ so that C x {t) = 
C2(t), Vi > 0. The condition cos (^3 = along with 
the conservation law ( |IL6[ ) then imply that smip(t) = 
0, Vt > 0. For si mplici ty we take ip(t) = 0, Vi > 0. 

Now, equati on (IIL5 ) determines R « \J\/2 — > 
but equation ( |IIL6 ) gives a divergent shift: A « 
log [4Ci(0) 2 /| J|J — > 00. However a consistent limit is ob- 
tained for the p hysical quantities. For example, from 



equations (III. 7 1 we obtain the amplitudes C\(t) 2 rj 
C 2 (t) 2 « (7i70)^exp(2Zr) . The evolution equation for 
r is still of the form 



dr 
dt 



E(0) - V lh 



(III. 18) 



lim V(t), which 

7^0 



but with a new potential Vn m (T) 
reads 

Vii ro (r) = Ci(0) 2 exp(2Zr) - 2Ft sin^ 3 . (111.19) 

Recall that cos <ps = . We conclude that the new 
potential Vy lm (r) is attractive if and only if sin </? 3 = 
sgn(ZF). This gives two separate cases, depending on 
the initial conditions: a bounded case and an unbounded 



1. The bounded sub-case of J = 0, cos^ = 



where we take k = for the minimum energy and k 
for the maximum. 



2. The unbounded sub-case of J = 0, costps = 

The unbounded sub-case is obtained if we set sin ip^ = 
— sgn (Z F) . In this case, the potential V[i m (T) in equa- 
tion (III. 19 1 has a slope of fixed sign. Depending on 



the initial conditions, the motion can have one turning 
point but then the variable Z r(t) will tend to —00 in 
the asymptotic way Z r(t) k, —\ZF\t 2 , so that the cor- 
responding amplitudes' squares will behave as C^{t) 2 « 
4F 2 t 2 , d(t) 2 = C 2 {t) 2 w Ke-\ zp \ t2 , corresponding 
to unbounded growth of the forced mode and fast decay 
of the other two modes. 



IV. PERIODS AND MAXIMUM AMPLITUDES: 
PARAMETRIC STUDY 

Having established for H = the general integrabil- 
ity of Eqs. (II. 4 1 and the boundedness of the solution for 



almost all initial conditions (with the exception of one 
sub-case), we turn our attention to two important phys- 
ical questions: 

(i) How does the physical period of oscillation depend 
on the parameters Z, F and on the initial condi- 
tions J, (fs, E(0) of the problem? 



The bounded sub-case is obtained if we set sin ^3 
sgn(ZF). Although there is no simple solution of the 



system (III. 18 1 in terms of elliptic functions, it is illustra- 



tive to obtain explicit expressions for the turning points, 
and consequently for the maximum and minimum values 
of the total energy (111.1), in terms of product log func- 
tions. The equations Viim(T m i n ) = Vu ro (r max ) = E(0) 
have solutions 



min / max 



^(O)sgn(^) 
2\F\ 

1 / |Z|d(0) 2 

— Wh — J — — r^- exp 
2Z k \ F v 



\F\ 



where k £ 1 denotes the different branches of the prod- 
uct log function Wf-(z), also known as the Lambert W- 
function. Notice that the argument of Wk appearing 
above is strictly negative in this case, restricted to the 
interval [— 1/e, 0] , with the lower bound being obtained 
only when 6*3(0) = 0. Since we require real solutions, we 
are limited in choice to the two branches Wo and W—\. 

Substituting these expressions for the turning points 
into equation ( III. 4 1 , wc find the maximum and minimum 



values of the total energy of the system: 



E 



min / max 



\F\ 



exp 



E(0)\Z\ 
\F\ 



(111.20) 



(ii) How do the upper and lower bounds on the total 
energy of the system behave as a function of the 
parameters Z, F and the initial conditions? 



A. Summary of analytic results 

Before embarking on a numerical parametric study 
of the behaviour of the physical period and maximum 
energy, it is useful to try to answer the above questions 
in the few sub-cases where an explicit solution has been 
found: 

• In the case H = 0, sin<p 3 = and J 7^ 0, as detailed 
in Section [ill C| the full solutions for the amplitudes and 
phases are available explicitly. In particular, the period 
of oscillations is given by 



T : 



2K(m) 



Zy/E(0)+R' 



where we have defined 




z 2 



and K(m) is the complete elliptic integral of the first 
kind. 

In this case the total energy ( |III.l I of the system 
happens to be constant, equal to E(0). This can be 
clearly identified in Figs. ([2| & ^ where the surfaces 
for maximum and minimum energy coincide on the line 
given by sin ip 3 = 0. 

• In the bounded sub-case of H = 0, cos ips = 
and J = 0, as detailed in Section |III D 1| although we 
do not have explicit expressions for the solutions for the 
amplitudes or even for the period of oscillations, we do 
have an exact expression for the maximum and minimum 
values of the total energy, from equation ( III.20[ ): 



mm / max 



w\ Wk l — ^T~ exp 



E(0)\Z\ 



where we take k = for the minimum energy and 
k = — 1 for the maximum, and we have assumed for 
simplicity that Ci(0) = C*2(0) and ip = 0. 

We now ask how these bounds obtained for the energy 
behave in the limit of arbitrarily large forcing. Recall 
that the Hamiltonian in the case cos ^3 = reads 

H = 2ZC 1 C 2 C 3 smtp. 

Assuming each amplitude to be non-zero, fixing H = 
means that we have had to set tp = 0, irrespective of 
the value of F. Therefore, by choosing appropriate initial 
amplitudes/phases so that J = 0, (f3 = ±ir/2 and H = 0, 
we can consider directly the limiting case when F tends 
to infinity but with fixed initial conditions. Since Wq is 
differentiable at zero, we have 

lim £ min = d(0) 2 . 

F-v+oo 

Calculating an equivalent expression for £ max , is not 
quite so straightforward since W-\ exhibits a singular- 
ity at 0. However, knowing the asymptotic behaviour of 
W-i(x) as x — > Q~ [26] , we find that to leading order as 
F +oo, 



E n 



\F\ 



In 



\F\ 



|Z|d(Q) s 



£7(0). 



(IV.l) 
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FIG. 2. -Emax and -E m i n as a function of </?3 and J. As a choice 
of parameters and initial conditions, we have fixed Z = 1, 
E(0) = 3, C 3 (0) = 1 and F = 1. Values of Ci(0), C 2 (0) and 
ip(0) can be found directly from the two conservation laws J 
and H = 0. 



in the maximal energy in response to forcing, in a similar 
fashion as shown in the asymptotic result (IV.l I. 



Period of oscillations. Since we only have an explicit 
expression for the period when sin (p 3 = 0, this too must 
be calculated numerically. After calculating numerically 
r max and T m ; n , integrating Eq. (III. 9 1 gives the period 



T = 2 



y/E(p)-V(u) 



This integral can be evaluated numerically, with special 
care being taken where the integrand has singularities 
at the two end-points. The results are shown in figures 
HandEl 

One thing to note with regards to the figures, is that 
setting H = restricts the possible range of parameters 
that can be taken. The effect of this can be seen in 
figures [3] and [5j where forbidden parameter values are 
indicated by 'hashed' regions at the base of the plot. 



B. Numerical study for arbitrary values of J and ^3 



APPROXIMATE SOLUTIONS IN THE 
INTEGRABLE CASE H = 



Maximum energy. When neither of the two analytic 
cases detailed in the above bullet points apply, the two 
roots r m i n and T max must be computed numerically from 
equations 15 ) . The resulting values of £7 max against 
varying ip 3 and J under fixed forcing are shown in 
figure [2j The unbounded sub-case as J — > 0, cos </? 3 = 
can clearly be identified. If instead we fix J > and 
allow the strength of the forcing F to vary, we obtain 
figure [3] Here we observe a less-than-quadratic growth 



As we have seen in Sec. |III| in the integrable case H = 
only certain choices of parameter values lead to explicit 



formulae for the solution of equation (III.9), in terms of 
elliptic functions. For other parameter values it is un- 
likely that such formulae exist and we have to resort to 
finding an approximate solution for r(t). Since the poten- 
tial V(t) = i?cosh (2Zt + A) — 2Ft sin ip^ is attractive, 
it is possible to approximate this potential by a polyno- 
mial in t. A quartic polynomial is the best compromise 
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FIG. 3. -Bmax and -E m i n as a function of ips and F. As a choice 
of parameters and initial conditions, we have fixed Z = 1, 
E(0) = 3, C 3 (0) = 1 and J = 0.1. Values of d(0), C 2 (0) and 
tp(0) can be found directly from the two conservation laws J 
and H = 0. 



FIG. 5. Period as a function of ipz and F. As a choice 
of parameters and initial conditions, we have fixed Z = 1, 
£(0) = 3, C 3 (0) = 1 and J = 0.1. Values of C*i(0), C 2 (0) and 
ip(0) can be found directly from the two conservation laws J 
and H = 0. 




FIG. 4. Period as a function of <^3 and J. As a choice 
of parameters and initial conditions, we have fixed Z = 1, 
£(0) = 3, C 3 (0) = 1 and F = 1. Values of d(0), C 2 (0) and 
tp(0) can be found directly from the two conservation laws J 
and H = 0. 



between simplicity and accuracy: simplicity in that the 
resulting approximate evolution equation for r(t) can be 
explicitly integrated, and accuracy in that one can cap- 
ture nonlinear effects such as the dependence of the pe- 
riod on the parameters of the problem. 

The problem is then reduced to compute the coeffi- 
cients of this approximate potential. A natural place to 
start would be to obtain these coefficients directl y from 
the fourth-order Taylor expansion of equation (111.9 1, 
about the minimum of the potential r = r*. However, this 



'local' approach is doomed because it does not eliminate 
secular terms, and therefore the solution of the approx- 



imate equations differs more and more from the exact 
solution after each period. 

Instead, the coefficients of this approximate potential 
must be chosen by a procedure that takes into account 
the global character of the motion, i.e., the fact that the 
variable r(t) can oscillate with high amplitude. 

A quartic polynomial in r has five free coefficients. We 
will determine two of them by requesting the approxi- 
mated potential to have the same turning points as the 
exact potential. Thus, our approximate system is defined 
as follows: 



dr 
dt 



= (t - r min )(r - T max )(a 2 T 2 + air + a ) 



wher e r m j n , r max are the turning points, defined in equa- 
tion (III. 12 1, and ao,a\,a2, are some unknown coeffi- 
cients. This approximate system replaces the exact evo- 
lution equation (III. 9). 



The remaining three coefficients will be determined by 
a minimisation procedure of the L? norm of the difference 
between the exact potential and the approximation. 



A. Approximate system, its solution and 
comparison with numerical simulations of the exact 
system 

In order to determine the unknown coefficients 
uq, a±, OL2, we impose a procedure that minimises the L 2 
error, globally over the interval [r min , r max ]. Let us define 
the function to be minimised: 



$(a ,ai,a 2 ) = 



1 



(P(r) - E(0) + V{t)Y dr, 
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where P(r) = (r - r lnin )(r - r max )(a 2 r 2 + ajT + a ) and 
V(r) is the potential energy, defined in equation (III. 10 1. 



The function $ is a positive-definite quadratic function 
of (ao,ai,ct2) € R 3 . It has a single minimum, attained 
at the point (ao, ai,a%) — (ctQ,a*,a%) € K 3 . This point 
is the solution of the linear system 



<9$ 

3a, 



(a* ,al,a* 2 )=0, j = 0,1,2. 



More explicitly, defining the real coefficients 

7j = y -r J (T-r min ) 2 (T-r max ) 2 dr, j=0,...,4, 

<5j = J r J (T-r min )(r-r mllx )[£;(0)-y(T)] dr, .3=0,1,2 

(elementary functions of T m ; n and T max ), we obtain the 
system 

' 7o 7i 72 

71 72 73 

72 73 74 

Since the matrix consisting of all the 7-coefficients is 

the solution for a* is 



1 




(8 


at 


- 


Si 






\ s 2 



< T„ 



invertible as long as r„ 
uniquely determined. 

After calculating the coefficients a* for j = 0, 1,2 as 
outlined above, we return to the approximate nonlinear 
differential equation 



dr 
dt 



(r- 



x)(r- 



c )(a* 2 T 2 + atr + a*). (V.l) 



Solving this equation for r(i) is now a matter of cor- 
rectly evaluating an elliptic integral. To do this we 



must first reduce equation (V.l) to Legendre normal 
form. The method closely follows the approach detailed 
in [37] : the quartic polynomial above is written as a prod- 
uct of two quadratics, Qi(t) = a* 2 r 2 + a*T + aj and 
Q 2 {t) = (t - r min )(r - r max ). We look for values of A 
such that Qi — XQ 2 is a perfect square of a first-degree 
binomial in r. To that end we require that 



4 («2 - ^)(«0 _ AT lni „T max ) ~ («! + (r min + T max )A) 2 = 0. 

Solving this quadratic equation for A we get 



Ah 



-u)± ^lu 2 - ( 



Tm ax Trr 



((at) 2 - 4a* a* 2 ) 



( T n 



where we have defined 



2«q + a*(r max + r n 



2a2 T max T m j n . 



Having found A + , A_ we can then write 

a* + (T m ax + Tmin)A- 



Qi - A_Q 2 = ("a - A_) r + 



<3i - A+Q 2 = («a ~ A+) t 



2K-A_) 

*** + (Tmax + T min )A + 

2(a* 2 - A+) 



Solving this system for Qi and Q 2 we obtain 



A + (Q!2 - X-) , 2 

n = — i 5 (t + a) - 



A-l - A_ 



a^ — A_ 2 a 2 — /\+ , 

Q2 = ^ —{T + a) — (t + P) , 



2 _ X-(a* 2 - A + )^_ t ... 
A 

aS — Aj 



A, - A_ 



'(r + /?) 2 , 



A, - A_ 



A+- A_ 



where we have defined the two constants 

_ a l + ( r max + T min)A- 



2(a* 

a l "T" (T ma x 



a-; 

f Tn 



0A H 



2(a* 2 - A+) 

The method will work if and only if A± are real and 
different. Since the coefficients a*, u, r max , r min are all 
real, the condition for the method to work is equivalent 
to the condition us 2 - (r ma x - r min ) 2 ((a 1 ) 2 - 4apa2) > 0. 
In turn, this inequality can be restated in terms of the 
zeroes r + , t_ of Q\ and r max , r m ; n of Q 2 , as follows: 
(Tmax - t+ ) (Tmax -T_)(Tmi n -T + ) (rmin -t_) > 0. Graph- 
ically, the latter inequality is equivalent to the statement 
that the zeroes of Qi are either: (i) complex (which neces- 
sarily come in conjugate pairs), (ii) real, but the interval 
(TminiTmax) contains either both zeroes of Qi or none of 
them. 

Can we establish that the method will always work 
for physically admissible values of the parameters Z, 
F and initial conditions R, (p^, E(0)7 Let us assume 
that r+, t_, the zeroes of Qi, are real. By looking 
at the minimisation problem, notice that the potential 
V(t) is convex. Therefore, if any of the zeroes r + , t_ 
were in (T max , T m j n ), then the approximating polynomial 
P(t) = (r— T m in)(T— T ma x)Qi(T) would have a zero inside 
the interval (T max ,T m j n ), and so P(t) would not be con- 



cave on r. 



mm j ' maxj 



It follows that a deformation of the 
coefficients aj could be found so that $(ao, a±, a 2 ) would 
decrease, thus violating the minimum principle that gen- 
erated the point a*. 

In practice, in numerical applications and for a wide 
range of parameters we have found that the zeroes of Q\ 
are complex (and therefore, conjugate pairs). In this case 
the method simplifies substantially because the quanti- 
ties A± satisfy A_ < < A+. 

So let us assume from here on, and for simplicity, that 
A_ < < A + . If we now make the substitution r = 
(t + a)/(r + p), we can reduce equation (V.l) to the 
normal form 



dr 
dt 



X+(a* 2 



-(/?-a) 2 (a 2 + r 2 )(r 2 -6 2 ), 



(A+-A_) 2 

where a, b are positive parameters defined by 

(ga - A+) \-(a* 2 -\+) 
(a* 2 -X.)' A + (a*-A_) 
("2-A+) A_(aJ-A + ) 



b 2 = - s 



(a* 2 -X-Y A + (a*-A_) 
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For physically sensible solutions (i.e., r(t) £ M) we must 
choose the branch r(t) 2 > b 2 for all t > 0. We obtain, in 
terms of the Jacobi elliptic function, 



too 



r{t) = b nc f u; t + u 

defined by the physical period, 

Tq ^ 2K(m ) 

and elliptic modulus 



to 



a 

i 2 + b 2 ' 



The value of the shift uo can be found directly in terms 
of the incomplete elliptic integral of the first kind: 



uo = siF ( arccos 



b(A + /?) 
A + a 



m 



with si G { — 1, +1} being given by 

si=agn(Z<7 3 (0)G8-a)). 

Undoing the transformations applied to r we finally un- 
cover the solution 



r(i) 



bp — a cn 1 — f-^t + u ° 




I 4K(m ). , 

cn 1 — j^-t + u o 




-b 



(V.2) 



Approximate expressions for the squares of the ampli- 
tude s Gj(t ) 2 are obtained directly from equations (III. 7) 
and (III.8), and the dynamical phase is obtained from 



equation (III. 11 1. Plugging the solution (V.2) of our 
approximate system, into these expressions, we gener- 
ate plots that can be compared with the numerical so- 
lutions of the original system (II.7). The results, shown 
in figures 6(a) to |6(d)| show excellent agreement. Wc 
have validated numerically, for a wide range of choices 
of parameters, the following condition of boundedness of 
the approximate solution (V.2): b > 1 or, equivalently, 
a* 2 < 0. 



VI. POINCARE SECTIONS AND EVIDENCE 
FOR BOUNDEDNESS WHEN H 

When the Hamiltonian is non-zero, we have no a priori 
guarantee over integrability of the forced triad system, 
due to the system's formal equivalence to an autonomous 
volume-preserving 3-dimensional system of first order 
equations. If we try to reduce the forced triad system 
to an explicit 3-dimensional system we face yet another 
difficulty: the generic appearance of coordinate singular- 
ities, as exemplified in equations (II. 2) at C1C2 = 0. 





(a)C 3 (t) 2 with F = 1 




(b)C 3 (t) 2 with F = 10 



o ir\ o, r 



(c)tp(t) with F = 1 



(d)yj(t) with F = 10 



FIG. 6. Comparison of the L 2 -method approximate solutions 
for the square amplitudes and dynamical phase. The cor- 
responding numerically generated Runge-Kutta solutions are 
also included as a direct comparison against different forc- 
ing strengths. Choices for parameters: 953(0) = —0.99^, 



E(0) = 3, J 
H = 0. 



1, Z — 1 and <p(0) determined by setting 



Fortunately, we can solve this difficulty completely 
and simplify matters considerably through an intelligent 
choice of variables. The procedure works only when 
H 7^ 0, and is based on the observation that 953 is mono- 
tonic and so it becomes a proxy for time. Let us define 
the new variables, 

p = $t{B 1 B 2 }, q = <S{B 1 B 2 }, 

with tps — Arg B3 as before. Using the conservation law 
J, we observe that 



(1^1 



\Bn 



J 2 = A\B X 



l 2 |S,l 2 , 



and rewriting the right-hand term \Bi\ 2 \B2\ 2 = p 2 + q 2 , 
we derive the remarkable formula 

|£i| 2 + \B 2 \ 2 = ^JJ 2 + A{p 2 + q 2 ). 

Differentiating the new variables (p,q,y>3) with respect 
to time we obtain, after straightforward computations, 
the following three-dimensional system: 

= J 2 +A(p 2 +q 2 )\B z \ cos p 3 , 



P 

q = Z^J 2 ~ 
H 

^3 = — onn 



4(p 2 + g 2 )|B 3 |smp3, 



where the amplitude I-B3I is a simple function of (p, q, 993) 
and the Hamiltonian H, which now reads 



H = 2Z\B* 



cos v?3 — p sin ip3 



(VIA) 
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In this representation, we see that the angle ip 3 is not 
only monotonically increasing or decreasing depending 
on H, but can also be read directly as the slope of the 
parametric curve (p(t),q(t)). This follows simply from 
the relation 



dq 
dp 



tan i/33. 



Furthermore, we see from this new representation of H 
given in Eq.(VI.l), that the condition for \B^\ to re- 
main finite is equivalent to the statement that the cross- 
product between the tangent vector, (p, q) and the off- 
centre radius vector (p, q — F/Z) must remain non-zero. 
In essence, (p, q) x (p,q — F/Z) must point out of the 
plane. 

Using this equivalence it makes sense to shift the origin 
of (p, q) coordinates and define a new pair of variables 
(r, (j>) such that p — r cos <j> and q = F/Z + r sin </>, or 

r 2 = (q-F/Zf +p\ 
q-F/Z 



tan d> — 



P 



Differentiating these expressions with respect to time 
gives the following system: 

H 



— VJ 2 + 4(p 2 + q 2 ) COt (cj) - tp 3 ) . 



<^3 = 



H 

2^2 
2Z 2 r 2 



sin 2 (<f>-(f3) 



The Hamiltonian now reads 

H = 2Z\B 3 \rsm (<p-ip s ). 

We immediately see that cj> as well as ^3 must be mono- 
tonic, increasing or decreasing identically with ip 3 . Since 
we have defined r to be positive, we also see that since we 
must necessarily have I-B3I > then either sin (cf) — ^3) > 
when ZH > 0, or, sin (<j) - <p 3 ) < when ZH < 0. De- 
pending on the choice of initial conditions, we therefore 
obtain the bounds 

jo < 4>- tp 3 <tt, iiZH>Q- 1 
j-7T < tp 3 < 0, if ZH<0. 

These bounds cannot be saturated, otherwise we violate 
the assumption H =/= 0. The variable cj> is monotonic, 
so we will define it as our new 'time'. Introducing the 
variable = cj) — 993 we get the new system 

dr 

— = —rcotQ, 
d(j> 

d9 4Z 2 r 4 sin 2 6> 

— = 1 



which is non-autonomous, but one where the 'time'- 
dependence appears only in the J 2 + 4(p 2 + q 2 ) term. 
Writing this term in these latest variables reads 



J 2 +4(p 2 + q 2 ) = 

where we have defined the < 
c(4>) according to 



Fsini 



c{J>f 



-dependent, positive function 



4 



F 2 cos 2 



One-dimensional particle interpretation. One may 

ask what is gained by these successive transformations. 
After all, this newest version of the system appears to 
offer little improvement over the original one written in 
terms of the variables p and q. Surprisingly, however, 
this new system has a novel one-dimensional particle in- 
terpretation similar to that obtained in the H 7^ case 
described earlier. If we make one final transformation 
h{4>) — r(0) _1 , we derive the one-dimensional particle 
representation: 



d 2 h 

ckA 2 



h = 



2Z 2 



H 2 h 3 



(VI.2) 



F sin 
Z 



c#) 2 



This can be interpreted as a driven harmonic oscillator, 
with driving force that is 27r-periodic in <f> but also de- 
pendent on the position h. Alternatively, we can write 
the system in terms of a potential V(h, 4>): 

where the 'dot' now denotes derivative with respect to 
<fr, and subscripts denote partial derivatives. Integrating 
Eq. (VI.2) with respect to h, we evaluate this potential 




to be of the form, 



V(h,<f>) 



(VI.3) 

where f((f>) is a yet to be determined function that is 
dependent on <f> only. This potential is attractive for all 
values of 4>, from which it follows that the solution must 
at least be bounded for finite 'time'. 

The important question that remains is whether the so- 
lution remains bounded for all time. Our analytical work 
is in progress, based on Levi's theory [55], and will be re- 
ported elsewhere. Our numerical studies suggest that the 
answer is positive. We present now the evidence. Due 
to the periodic nature of the forcing term, we define a 
Poincare map, P : S — >• S, where S = M>o x M is the 
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Furthermore, since the map is continuous when restricted 
to each of these one-dimensional loops, Brouwer's fixed- 
point theorem guarantees the existence of at least one 
fixed-point. It therefore follows that the system contains 
infinitely many periodic orbits, with period 27r in <f), one 
of which can clearly be identified on the line = where 
these loops converge to a point. 



VII. CONCLUSIONS 



10 15 20 25 30 35 

m 



-4 , 



(a)H = 0.1 



( C 



IS 



m 

(b)H = 1 

FIG. 7. Poincare sections produced for varying values of H. 
Here we have fixed the values of Z = 1, J = 1 and F = 1 in 
both panes. Closed invariant curves can clearly be identified 
in each figure. 



set of initial conditions of the form (h(0), /ty(0)), and the 
transverse section is given by </>„ = 2mr corresponding to 
H < 0, and <f> n = -2nir if H > 0, taking n £ N . Suc- 
cessive iterations of the map P are then defined through 
the equation 

p n (h(o), mo)) = W^), M&0)> ™ e N o- 

Examples of some numerically generated Poincare sec- 
tions are given in figure [7] There is clearly no chaotic 
behaviour in the system, with the motion being quasi- 
periodic in nature, restricted to clearly visible closed in- 
variant curves. Each of these curves is in fact a loop, 
suggesting that the motion remains bounded indefinitely. 



We have presented an in-depth analysis of the triad 
equations driven by external forcing, Eqs. |II.1| We ob- 
tained a complete understanding of the dynamics when 
the Hamiltonian is zero by using conservation laws to 
reduce Eqs. to the one-dimensional motion of a par- 
ticle in a time-independent potential. In this case, the 
dynamics are integrable and bounded for almost all ini- 
tial conditions with a single exception provided by the 
special initial condition discussed in Sec. |IIID 2| Along 
the way, we presented several new explicit solutions, ana- 
lytic formulae for the period and maximum energy of the 
motion and a novel scheme for approximating the dynam- 
ics when an explicit solution seems beyond reach. Our 
results for the case H ^ are more qualitative and based 
on reducing the problem to the one-dimensional motion 
of a particle in a time-periodic potential. Poincare sec- 
tions of the dynamics provide strong empirical evidence 
that the amplitudes remain bounded when H ^ and 
are typically quasi-periodic although periodic orbits also 
exist. Broadly speaking, our results illustrate two im- 
portant points which are likely to have relevance beyond 
the particular problem studied in this article. The first 
is that, as illustrated by our boundedness results, the 
addition of forcing to a nonlinear wave system does not 
necessarily provide a continual source of energy for the 
system. This should provide pause for thought for re- 
searchers interested in performing numerical simulations 
of wave turbulence (where a constant source of energy 
is often required in order to compare with theory) since 
such simulations are often performed by adding an exter- 
nal forcing term to a Hamiltonian wave equation. The 
second general lesson to be drawn from this work is the 
importance of phases in understanding the dynamics of 
resonantly interacting nonlinear waves. Since the dynam- 
ics of phases are often more difficult to obtain theoret- 
ically, they are often somewhat under-emphasized. Our 
results present some of the most striking examples yet 
of how varying an initial phase can result in completely 
different behaviour. 
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